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The dynamical behavior of a kind of models with hierarchically constrained dynamics is inves- 
tigated. The models exhibit many properties resembling real structural glasses. In particular, we 
focus on the study of time-dependent temperature processes. In cooling processes, a phenomenon 
analogous to the laboratory glass transition appears. The residual properties are analytically eval- 
uated, and the concept of Active temperature is discussed on a physical base. The evolution of the 
system in heating processes is governed by the existence of a normal solution of the evolution equa- 
tions, which is approached by all the other solutions. This trend of the system is directly related 
to the glassy hysteresis effects shown by these systems. The existence of the normal solution is 
not restricted to the linear regime around equilibrium, but it is defined for any arbitrary, far from 
equilibrium, situation. 

PACS Numbers: 05.70.Ln, 64.70.Pf, 45.70.-n 



I. INTRODUCTION 



In recent years, there has been quite a considerable amount of work in models in which glassy behavior is generated 
not by quenched disorder, but by kinetic constraints. The kinetic restrictions are responsible for the slow relaxation, 
since the state of a particle or a group of particles can only change if some condition of its environment is fulfilled. In 
particular, "facilitated" models have been considered, both for structural glasses and for granular systems [[|. 
The characteristic feature of facilitated models is that a particle (spin) can only change its state if a certain number 
of its neighbors is in an excited state. Also, hierarchically constrained models have been used to study stretched 
exponential relaxation in glasses Q. In these models, the system is structured in levels and a particle in a given level 
can only make a transition if a given cluster of particles in the lower level is in a certain subset of configurations. 
Then, the dynamics of the several levels are coupled, and the characteristic relaxation times increase with the level 
index. Hierarchically constrained dynamics may be relevant for those complex systems in which the time evolution of 
the slowest modes is controlled by the relaxation of the fastest ones. This qualitative picture is adequate to describe, 
among other problems, proteins relaxation (§) and the densification of powders and structural glasses at high pressure 
P,^fj|. Very recently, a kind of hierarchically constrained dynamics has been shown to exhibit, in quite a natural way, 
logarithmic relaxation . This kind of "anomalous" , highly nonexponential, decay is observed in a wide variety of 
complex systems, including spin-glasses 12 granular materials P, p^| , p^ , p^[ , structural glasses Jl~6|— |l8|] and protein 
models |Q. 

The aim of this work is to study the dynamical behavior of a general class of hierarchically constrained models 
when submitted to more complicated processes. In particular, we are interested in the behavior of a system with 
hierarchical constraints when the temperature changes in time, what makes the coupling between the levels be time- 
dependent. The consideration of time-dependent temperature processes requires an extension of the original model 
as formulated by Palmer et al. 0. This will be done in a very simple, but natural, way: the coupling between the 
levels vary in time because the probability of a cluster configuration allowing to relax a particle in the next level 
depends on the temperature. On the other hand, neither the number of particles in a given level nor the length of 
the clusters "facilitating" the relaxation depend on the temperature. They are considered as quantities defined in the 
coarse-grained description of the system introduced to model the physical problem at hand. 

Let us note that hierarchical models can also be applied to the analysis of non-thermal systems, as granular materials 
in the dense regime. For those materials, thermal energy is not enough to make the system explore the phase space of 
configurations. Then, the system must be externally excited, for instance vibrating it, in order to be able to evolve. 
In these situations, the role of the temperature is played by the intensity of the external dri ving . If the stationary 
state reached by the system in the long time limit can be described by Edward's theory (2^J|l|, the compactivity 
X, which is the analogous of the temperature in thermal systems, will be a function of the intensity of the external 
force. Then, by exploiting the analogies of Edward's theory, i. e., substituting volume by energy and compactivity by 
temperature, it is possible to incorporate non-thermal systems in our formulation. 

Processes in which the temperature is time-dependent are physically relevant because they can be used to study 
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some characteristic dynamical aspects of glasses. For instance, when a supercooled liquid is cooled down to very low 
temperatures, a laboratory glass transition is observed. A dramatic change in the behavior of the system takes place, 
and it departs from the equilibrium curve, getting "frozen" in a far from equilibrium state. This transition appears 
as a consequence of the fast increase of the relaxation time with decreasing temperature. In order to characterize 
the cooling process, experimental physicists often use the residual value of the relevant physical properties, i. e., the 
difference between their actual values over the cooling curve and the value obtained by extrapolation of the equilibrium 
curve to the very low temperature region p^ , |23| . If the system is reheated from the nonequilibrium state, it returns to 
equilibrium for high temperatures, but it follows a different curve from the cooling one, giving rise to hysteresis effects. 
This phenomenon is related to the "nonlinearity" of glassy relaxation: the approach towards the equilibrium curve 
depends on the configuration of the system, measured by the so-called fictive temperature |22]-|2~4|| . Narayanaswami's 
theory provides a phenomenological explanation of this behavior [p2|-|25|] . Interestingly, a similar behavior has been 
found in vibrated granular materials when the tapping intensity is varied in a cyclic way |2q] , although the hysteresis 
effects are more evident when the heating process begins in a loosely packed state, referred as to the "irreversible" 
branch in the experiments. 

We will start from a very general hierarchical spin model, in which pseudospins are organized into levels, labelled 
by an index n. The pseudospins are assumed to correspond to some coarse-grained description of the system. They 
can take only two values, representing, for instance, two possible densities of a certain small subvolume of the system. 
One pseudospin in level n + 1 can only flip between the two possible values if a cluster of /i n spins in level n is in 
a given subset of configurations. This is the basic characteristic of hierarchically constrained models as introduced 
in Ref. 0, and it slows down the relaxation in level n + 1, as compared with that of level n. Here we will consider 
the simple choice that all the spins in the cluster must be in the up (excited) state. The exact dynamical equation 
for the evolution of the pseudospins involves very complicated moments of the probability distribution. To get an 
exactly solvable model, a "mean-field" approximation will be introduced. Then, the characteristic relaxation time r n 
of level n is seen to increase both with the index label n, due to the hierarchical constraint, and also with decreasing 
temperature, since the configurations allowing the system to relax become less probable when the temperature is 
lowered. 

Some exact dynamical results for systems described with master equations with time-dependent transition rates are 
known. In particular, the existence of a "normal" solution, i. e., a solution of the master equation that is approached 
by all the others, has been proved on a very general basis [ p7[ . The main required conditions are the irreducibility of 
the Markov process for long enough times and that the transition rates be externally controlled, so that they do not 
depend on the probability distribution of the system. Moreover, it has been established that the normal solution tends 
to the equilibrium curve for very high temperatures, in continuous heating processes |g7[| . The linear correction of the 
normal solution with respect to the equilibrium curve has been computed using Hubert's method |28|. It is not evident 
whether these results still hold when approximations are introduced in the dynamics of the system. In particular, in 
"mean-field" type approximations, the demonstrations in Refs. (27jj2^] are not valid, since the transition rates become 
functionals of the probability distribution function. Nevertheless, we will show that all the above mentioned properties 
apply to our simplified model. This is a good test of the plausibility of the approximations carried out, and perhaps 
an indication that the results obtained here are more general than the derivations in Refs. [p7p8| . 

The organization of the paper is as follows. In Sec. || the hierarchical model is introduced, and the exact evolution 
equation for the average spin is obtained. By introducing a mean-field approximation, this equation can be closed. 
Afterwards, a specific, but quite general, choice for the functions defining the model is made. This allows to introduce a 
continuous limit in which the relaxation of the system at constant temperature is solved in Sec. III. Time-dependent 
temperature processes are considered in Sec. IV, where the general solution for the evolution of the probability 
distribution and the average energy are obtained. The gener al sol ution is similar to the expression proposed by 
Narayanaswami on a phenomenological basis |2^,|24|,^5| . Section IV A is devoted to the analysis of Hilbert's expansion, 
which is valid in the very high temperatures regime. Cooling processes are addressed in Sec. [y] where, for the sake 
of simplicity, a concrete cooling law is studied, for which the r esidu al properties are analytically calculated. A 
qualitative analysis of the glass-like transition is presented in Sec. V A . It allows to give very good estimates of the 
residual properties, and leads to the introduction of the concept of fictive temperature in a very natural way. The 
behavior of the system when it is reheated from low temperatures is considered in Sec. VI. The main role played by 
the normal solution for the understanding of the hysteresis effects shows up. Moreover, the analysis clearly indicates 
that the relevance of the normal solution is not restricted to near equilibrium situations. Finally, a discussion of the 



main points in this work is given in Sec. VH 
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II. DYNAMICS OF HIERARCHICALLY CONSTRAINED MODELS 



In this section a general kind of spin models with hierarchically constrained dynamics will be introduced. We will 
focus on the evolution of the average value of the spin, which is supposed to be the relevant variable. For instance, 
in a thermal system it will be directly related to the mean energy. Then, let us consider a system whose degrees of 
freedom can be classified into levels, labelled by an index n = 0, 1, 2, . . . , n max . The degrees of freedom in level n will 

(n) 

be represented by N n pseudospins, a\ = ±1, i — 1,2,..., N n . The hamiltonian of the system is assumed to have 
the form 



H = ^E m !"' -W = i+^. (2.1) 

=0 i=l 



Note that mf^ is the occupation number of the "up" (+1) state of the corresponding site. In order to write Eq. 
( |2.lD we have supposed that there is no interaction between the pseudospins, but there is an "external field" h. For a 
thermal system, Ti. gives the energy of a given microstate of the system, while for a non-thermal system, like a powder, 



it could be interpreted as the volume of a given, mechanically stable, configuration of "grains" 20 23]. Using the 
terminology for thermal systems, the average value of the dimensionless energy per spin over the ensemble of systems 
considered is 



(H) _ 1 V^V^fn) 



n=0 i=l 



where 



^ = (mf)) = l±^l (2.3) 

is the probability that the i-th spin of level be in the up state, N — X)ri=cT ^ n 18 ^ ne total number of pseudospins, 
and the angular brackets denote statistical ensemble average. 

The system is considered to be in contact with a heat bath at temperature T, so that the equilibrium average value 
of the pseudospins does not depend either on i nor on n, and it is given by 

(*) e = (vf\ = -t a x&Q-) . (2.4) 

Here T* is a dimensionless temperature, T* = 2ksT /h, fc^ being Boltzmann's constant. For the sake of concision, we 
will drop the asterisk in the following. The above average value of the spin follows from the equilibrium probability 
for the "up" state of any spin, 

Pe=p£ )= ff X ' (2-5) 

In the limit of infinite te mpe rature, or zero external field, both states of the pseudospins are equiprobable, p e = 1/2, 



and (a) e = 0. From Eq. (2.2), it follows that the equilibrium value of s is 

e e =Pe- (2.6) 

For granular materials, the role of the temperature T is played by the compactivity (2^j2l[], which is linked to the 
intensity of the perturbation allowing the system to explore the configuration space. 

The dynamics of the model is formulated by means of a master equation with single-spin flip Glauber transition 
rates p9| . Let us consider the flip of a given spin . This transition connects a given configuration er of the whole 
system with the configuration r!- (t, where f6 is the operator which rotates the spin cj n \ keeping all the other 

(n) 

spins the same. The transition rate for the flip of the spin a\ in configuration er is 



o\- tanh 



T 



(2.7) 
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The characteristic relaxation rate a ™ of the spin cr> nJ depends on the configuration er of the system through the 
hierarchical constraint 

fei+A 1 n-l-l 

II *«rC»-i\ +1 , ( 2 -8) 

where 5« is the Kronecker delta. This expression implies that the spin needs, in order to flip, that all the spins 
belonging to a cluster of length /x n -i starting at a given spin fej of level n — 1 must be in the up (+1) state. Besides, 
the characteristic flip rate of the spin is the average of the characteristic flip rates of the spins belonging to the cluster 
determining its possibility of change. This restricting condition is schematically depicted in Fig. [I]. Note that the 
possibility of a given spin in level n to flip is restricted by the state of a set of clusters in all levels n' < n, the number 
of clusters involved in each level increasing as n' decreases. The hierarchical constraint implies that the configuration 
with all the spins in the down (-1) state is completely absorbent, i. e., the system does not evolve in time from that 
configuration. 



(n) / _\ 



(n-l) 



n-2 



n-l 



n 



n+1 



(n) 

FIG. 1. The framed spin in level n, a\ , has a nonvanishing probability of changing its state only if the framed cluster of 



spins in level n 
and fj,„-i = 5. 



1 are in the state shown in the figure, i. e. all of them up (+1). In this example, we have taken ki = i — 2 



We are interested in the time evolution of the average spin a 

»\ _ 



(n) 



which is given in Glauber dynamics by [f29| 



(2.9) 



and substitution of Eqs. (2.7) and (2 
d 



dt 



(n— 1) i \ in— l) / \ 

a kt [cr)al i+1 (a) 



(n-l) 



into this expression yields 

l/Mn- 



'(^-w.) n 



-l) 



(2.10) 



-1 for all i, n. This equation is rather involved, since it couples 



where we have taken into account that (a, 

the evolution of (cr l - rl ' ) ) to moments of the probability distribution containing an increasing number of spins of all the 

levels n' such that < n 1 < n. The levels < n' < n — 2 enter into the equation through the rates otj ^(<r). Then, 
we introduce at this stage a sort of "mean field" approximation for the transition rates, which is closely related, in 
spirit, to the seminal work of Palmer et al. Q]]. Within this coarsening dynamics, it is reasonable to expect that spins 
in level n evolve over a time scale quite larger than that characteristic of level n — 1. This means that spins in level 
ri — l change many times their state before a transition in level n takes place. Thus, we replace the product of the 

-l)^(n-l) =+1,... 



Kronecker deltas in Eq. (2 
that all the /i. 



by its average value, i. e., by the probability P( n ^ (cr^ n 



i spins of the given cluster be in the up state, 



» 



(n-l) 
p(n-l 



/ \ (n — 1 



{*)■ 



(n-l) 

V 



(n-l) _ 

ki " 



T ("-l) 



-1) 



-1=+1) 



(2.11) 



4 



Moreover, we will restrict ourselves to situations where there is spatial homogeneity within each of the levels, so that 
the dependence on the specific site considered in a given level can be dropped, obtaining 



a (™)(«x) = a {n - 1 \ ( r)P^- 1 \a { r 1) = +1, 
Iteration of the above relation gives 



(2.12) 



aW(<r)=a (0) ft P °Vi C 

3=0 



0) 



-1,' 



-1), 



(2.13) 



with being a constant, that characterizes the relaxation rate of the spins belonging to level n — 0, whose dynamics 
is not constrained. As determines the basic time scale, which is arbitrary, we will take = 1 in t he fol lowing. 
In the mean field approximation just introduced, the time evolution of the average value of the spin, Eq. ( |2~T0| ), takes 
the form 



d_ 
~dt 



n-l 

- n m 

3=0 



+1,- 



r«) = +1) (V»>) - (a) e ) 



(2.14) 



Now, a new approximation will be made. The probability P«) in Eq. (2.14) will be substituted by its equilibrium 
value. This would be exact in linear response around equilibrium, but it will be taken here as an approximation 
leading to the basic equation of our hierarchically constrained model, namely 



A n—l 



3=0 



(2.15) 



where p e is the equilibrium probability of any spin being in the up state, given by Eq. ( |2.5| ). Again, on physical 
grounds, this is a sensible approximation due to the separation of the characteristic time scales of the different levels. 
Because of the hierarchically constrained dynamics, spins in level n — l reach equilibrium over a time scale in which 
spins in level n have not begun to evolve. However, as a consequence of the last approximation, the configuration 
with all the s pins in the down state is no longer absorbent. 

From Eq. ( |2.15 ), an equivalent equation can be written for the evolution of the probability of the up state in 
level n, defined in Eq. (|2.3p, i. e., 



(n) _ 



dt P 

being a n the characteristic relaxation rate of level n, 



= P 9 e n , 



9n 



n-l 
3=0 



(2.16) 



(2.17) 



Equation ( j2.16|) implies that, due to the hierarchical constrai nts, t he spin relaxation slows down with increasing level 
n, since a n is a decreasing function of n, because p e < 1. Eq. (2.16) is the main result in this section. In the following, 
we will explore its implications, considering first processes at constant temperature in the next section, and cooling 
and heating processes in the remainder of the paper. 



III. RELAXATION AT CONSTANT TEMPERATURE 



For the case of constant temperature T, and therefore constant a n , Eq. (2.16) is easily solved 



P (n) (0)-p ( 



(3.1) 



Then, each spin relaxes exponentially to equilibrium with the rate characteristic of its level. 

For the hom ogen ous situations within each level we are considering, the dimensionless mean energy per spin |3C 
defined in Eq. (2.2) simplifies to 



5 



e(t) 



E 



(3.2) 



where w n = N n /N is the fraction of spins in level n, verifying X)"=cT Wn = 1- 
Use of Eq. (O) into Eq. (fO) yields 



7 l u 



n=Q 



e(t)=e e (T)+ >>„ P (n) (0)-Pe e"""', (3.3) 



for the relaxation of the energy at constant temperature. In order to proceed, we will consider the simple case in 
which the initial probability distribution p^ (0) does not depend on the index level n. This will be the situation, for 
instance, when the initial state corresponds to equilibrium at a different temperature T + AT. Thus, the relaxation 
function of the physical property described by the hamiltonian of the system is given by 

v ' n— 

This equation is a generalization of the result derived by Palmer et al. in their pioneering work in hierarchically 
constra ined dynamics B , which corresponds to the choice p e — 1/2. This is formally equivalent to the particularization 
of Eq. (3.4) for T — > oo. If other positive values of the temperature are considered, the effect is an increase of the 
relaxation times 

r n = a; 1 , (3.5) 

since p e is a decreasing function of the temperature. A main advantage of the formulation of the hierarchical models 
as presented here, aside from its larger generality, is that it allows the analysis of processes in which the temperature 
of a thermal system (or the vibration intensity in a granular system) changes in time. This kind of processes will be 
addressed in the next section. A mean relaxation time r can be defined as 



dt<j)(t) = ^w n r n , (3.6) 



n=0 

providing a quantitative measure of the time it takes the system to relax to equilibrium at temperature T. 

Let us consider that the fraction of spins in level n, w n , and the number of "facilitating" spins in level n, fi n , depend 
very smoothly on n, i. e., they can be expressed as functions of the form 

w n = w{nrj) [i n = ii{nrj) , (3.7) 

where r\ <C 1. These seem to be sensible conditions when modelling a real system, in which the introduction of the 
levels and the pseudospins is associated to some coarse-grained description. By defining 

x„ = m) , (3.8) 

that is a continuous v ariab le in the limit r\ — > 0, the sums over n can be replaced by integrals. The relaxation rate of 
level n, given by Eq. ( 2.17 ), becomes a function of the continuous variable x, 

a{x) = pfW , (3.9) 

with 

Ek=o w kV J dx'w(x') 

where x max — tw i), and the normalization of the weights w n has been used. Therefore, the relaxation function <j)(t), 
given by Eq. (3.4), becomes 

m - - ' , (,n> 

Ln=o VW* Jo dxw(x) 



G 



in the continuous limit. 



In general, Eq. ( 3.11 ) is mathematically rather involved, since it depends both on the functions n(x) and w(x). The 
simplest possibility appears to be n(x) proportional to w(x), i. e., \x n proportional to w n or, equivalently, to N n . In 
other words, the simplest kind of hierarchically constrained models shows up when the number of "facilitating" spins 
at a level is an extensive function of the number of spins at the same level jjT| . This condition is expressed as 



/j,(x) = a w(x) , 

with a being a constant, independent of x. In this case, it is useful to define the new variable 

f„ x dx'w(x') 

u = ±M - - 

/ Cmax dx'w{x') ' 



(3.12) 



(3.13) 



measuring the fraction of spins belonging to levels up to n — x/rj. In terms of m, the relaxation rates of Eq. (3.9) are 
given by 



a{u) = p g e {u) , g(u) = au, 



the relaxation function is expressed as 



and the mean relaxation time reads 



dt<f>(t) = / duriu) = ^f- 

Jo a|lnp e 



(3.14) 



(3.15) 



(3.16) 



with t{u) — a -1 (it). It is interesting to consider situations for which p°\ <C 1, so that the minimum relaxation rate 
a(l) is much smaller than the maximum one a(0) = 1 in our dimensionless time scale. In this case, the relaxation 
function <f>{t) is linear in Int over an intermediate time window, 1 ^ i P7 a , namely [filj 



i - 



i 



a\ \np e 



(7 + In*) 



(3.17) 



where 7 stands for Euler's constant, 7 ~ 0.577. This kind of linear logarithmic behavior is characteristic of a great 
variety of complex systems, including spin-glasses Jl2],[l3]], granular materials j^,|ll],|l4|,|l5| , structural glasses [jig-fill , 
and protein models jpj[L9|. In the present conte xt, th e condition p° 1 corresponds to a "low" temperature limit, in 
which the mean relaxation time, given by Eq. (3.16), is very large, i. e., the relaxation of the system becomes very 
slow. 

It is worth noting that, in terms of the u variable, the continuous limit is formally obtained by changing the functions 
of the index level n by the corresponding functions of the variable u, and by making the replacement 



n=0 



du . 



(3.18) 



It follows that, in the continuous limit, the dynamical behavior of the system does not depend explicitly on the level 
populations N n , but only on the relaxation rates a expressed as functions of u, as given by Eq. ( p. 14 ). 



IV. TIME-DEPENDENT TEMPERATURE PROCESSES 

In this section processes in which the temperature changes in time will be studied. The evolution equation ( 2.16| ) 



is now 



j t pM(t) = -a n (T)[pW(t)-p e (T) 



(4.1) 



where T = T(t). The general solution of Eq. ([hi]) is 
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pW(t)= p(")(t )-Pe(To) X n (*,«0)+Pe(T)- / ttt 

to 



.d Pe (T')dT' 



dT' dt 



T Xn{t,f) 



(4.2) 



Here T = T(t ) is the initial value of the temperature, T = T(t), T = T(t'), T" = T(t"), and we have introduced 
the function 



- f* 1 dtct n (t) 



(4.3) 



The above equation is valid for any law of variation for the temperature. Taking into account Eq. (|j), the average 
energy per spin e, is given by pof 



n max 

:(t) =Y,w n [p(")(to) -Pe{To)] Xn(*,*o) 



-£e(T) 



t.de e (T')dT 



dt 



to 



dT' dt 



(4.4) 



n=0 



where e e is the average equilibrium energy, defined in Eq. (2.6). 

Let us assume that the system is initially at equilibrium with T = Tq. Then the first term on the rhs of Eq. (4.4) 
vanishes, and 



e{t)=e e {T) 



dt 



,de e {T') dT' 



to 



dT' dt 



T M{t,t'), 



(4.5) 



where 



M(t,t') = £ W nX n(t,t'). 



(4.6) 



n=0 



is a memory function. In the case of con stan t temperature, M(t,t') is e qual to the relaxation function cf>(t — t'), 
as seen by comparing Eq. ( |4.6| ) with En. (3.4). The structure of Eq. (4.5) is the same as that of Narayanaswami's 
phenomenological theory of glasses []22[- 25[ . A similar result was obtained some years ago for the one dimensional 
Ising model with Glauber dynamics PjJ. 



A. High temperature limit: Hilbert's method 



We are going to look for a solution of Eq. (4.1) by means of Hilbert's method. A special solution 



(n) 



Ph 



k=0 



is constructed in an iterative way as follows. We take 

pP>°(t)=p e (T), 

while for k > 1 



p { ^\t) = - a -\T) 



dt 



(4.7) 



(4.8) 



(4.9) 



Equation (4.8) shows that Hilb ert's expansion agrees with the equilibrium distribution to the lowest order. Besides, 
for k = 1 we get from Eq. (4.9) 



PP' 1 ® = ~r n (T) 



dp e (T) dT 
dT dt 



(4.10) 



This equation indicates the main limitation of Hilb ert's method. Due to the divergence of the relaxation times 
t„ = a^ 1 in the low temperature limit, see Eq. ( 2.17| ), also pW' 1 diverges in that limit. As a consequence, Hilbert's 
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solution is only accurate in the high temperature regime, in which an expansion around equilibrium provides a good 
approximation. Restricting ourselves to high temperatures, we approximate 



pW(T)~p B (T)-r B (T) 



d Pe (T) dT 
dT dt 



and, from Eq. (3.2), 



e H (T) ~ e e (T) 



de e (T) dT \ 



dT dt 



n=0 



W„T n (T) . 



(4.11) 



(4.12) 



Taking into account the definition of the average relaxation time t, Eq. (3.6), the above expression is seen to be 
equivalent to 



e H (T) ~ e e (T) 



de e (T) dT r 
dT ~dt' 



CO, 



(4.13) 



which agrees with the high temperature behavior of Eq. (4.5). In Eqs. (4.11) and (4.13), p^ and eh depend on time 
only through the temperature T(i). Thus, Hubert's method provides a "normal" solution, in the sense often used 



in kinetic theory. The validity of a expression identical to Eq. (4.13), also in the high temperature limit, has been 
established for a quite general class of systems whose dynamics is described by a master equation ^8|. Although 
we have made here several drastic approximations in order to get a closed equation for the average spin, the high 
temperature limit of the solution, given by Hilbert's method, remains formally the same as that of the exact solution 
of the original model. Certainly, this is a good property of those approximations. 

What is the physical meaning of the failure of Hilbert's expansion for low temperatures? Due to the divergence of 
the characteristic relaxation times, the system does not have enough time to relax to the equilibrium curve at very 
low temperatures, and it gets "frozen" in a far from equilibrium state. Since Hilbe rt's method is an expansion around 
equilibrium, it fails in the low temperature region. In fact, the rhs of Eq. (4.11) becomes negative for low enough 
temperatures. On the other hand, Hilbert's expansion is useful to estimate the values of the physical properties in the 
"frozen" state |28 . Also, Hilbert's method provides a qualitative understanding of the hysteresis effects appearing in 
thermal cycles (cooling and reheating). In cooling processes (dT/dt < 0), it is Eh > £e, while in heating processes 
(dT/dt > 0) it is Eh > E e . Then, eh lies to opposite sides of the equilibrium curve for cooling and heating processes, 



and hysteresis effects show up in thermal cycling experiments, as it will be discussed in more detail in Sec. VI 



V. COOLING PROCESSES 



Next, we are going to study the continuous cooling of the system down to very low temperatures. The origin of 
time is taken at the beginning of the cooling process. The initial condition will be the equilibrium configuration at a 
"high" temperature To, i. e., 



pW(0)=}>e(To). 

Then, the first correction in Hilbert's expansion, p^ 71 '' 1 ^), is very small as compared with p e (To) for T 
Particularization of Eq. (4.2) for the above initial condition gives 



dT' dt' 

Since p e (T) is an increasing function of the temperature, for continuous cooling processes it is 

p [n \t) > p e (T) for all t and n. 



(5.1) 
- T . 

(5.2) 
(5.3) 



The possible deviations from the equilibri um distribution always lead to an increase of the probability of the spin 
being in the excited state. Moreover, Eq. (|5.2|) directly implies that 



e{t) = e e (T) - J* dt'^l^Y. w " Xn(t, > e e (T) . 



(5.4) 



n=0 
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where we have used Eq. ( |3.2| ) . This inequality has been experimentally observed in glass- forming liquids [|22|,|23| . Since 
the reported experiments were made at constant pressure, the quantity e considered here must be interpreted as the 
enthalpy in that context. 

In order to proceed further in our analysis, the continuous limit introduced in the study of the relaxation at constant 
temperature in Sec. [II will be considered. As already mentioned, this continuous limit is expected to be closer to 
the description of real systems than the discr ete level picture. Besides, for the sake of concreteness, we will restrict 
ourselv es to those models verifying Eq. ( |3.12 ). The index level n is substituted by the continuous variable u, defined 
in E q. ([3.13 ), representing the fraction of the total number of spins up to level n. With an obvious change of notation, 
Eq. (5.2) becomes 



p(t;u) =p e (T) 



(i 



* , ,dpJT') dT' . . 

dt' - v ; x(t,t ;u) . 

dT' dt 1 AV ' ' ; 



where 



x(M>)^>'' a(TV \ 



with a(T;u) given by Eq. ( 3.14 ), i. e., 



a(T;u)= Pe (T)° 



Also, using Eq. (3.18), it is found 



e(t) = e e (T) 



dt 



,d e e (T') dV 
dT' dt' 



du x(t, t ;u) . 



(5.5) 

(5.6) 
(5.7) 

(5.8) 



The time evolution of the probability p(t; u) depends on the explicit form of the cooling law. In experiments, linear 
cooling is usually employed, 



dT 
dt 



(5.9) 



where r c > is the cooling rate determining the time scale, r c , over whi ch t he te mper ature changes. Linear cooling 
implies that p e (T) depends on time in a rather involved way. From Eqs. (2.5) and (5.9) one gets 



dp e 
dt 



-r c p e (l -p e ) In 



1-Pc 



(5.10) 



We are interested in cooling processes for which the temperature changes slowly in time, r c < 1, so that the system 
departs from the equilibrium curve for very low temperatures, where p e <C 1. Then, a law equivalent to linear cooling, 
aside from logarithmic corrections, is 



dp e 
dt 



-r c Pe , 



(5.11) 



having the advantage that analytical calculations are much more simple with Eq. ( 5.12 ) [|3l|-|33|. In the following, we 
will use the notation 



Pe(T)=p e , p e (T')=p / e , p e (T )=p e0 . 



(5.12) 



The time integrals in Eq. (5.5) can be transformed into integrals over p e by means of the cooling law ( 5.11 ), with the 
result 



p(t;u) =p e + 



dp' e exp 



,/aii au 

Pe Pe 

r c au 



(5.13) 



The second term in this expression is dominant in the low temperature region, where p e is very small. Therefore, 
it follows that spins in any level u fall in a nonequilibrium state for low enough temperatures. The details of this 
glass-like transition will be analyzed below, in a separate subsection. 

One of the main quantities characterizing a given cooling process is the residual value / ros of a relevant property /. 
The residual value measures the excess with respect to the equilibrium curve, extrapolated to very low temperatures, 
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/res = lim (/ - f e 



In particular, for the probability p(t; u) we have from Eq. ( 5.15 ) 



rPeO / plan 
Pres(w) = / ^exp — 

Jo \ r c au 



which is easily transformed into 



1 i r x o 

p res (u) = — {r c au) a 



dxx^ 1 e x . 



(5.14) 



(5.15) 



(5.16) 



where xq = p^q / {r c au) . The slow cooling limit is defined by the residual properties being independent of the initial 
conditions, and determined univocally by the cooling rate |31-34|]. In our case, slow cooling means that the upper 
integration limit in Eq. (5.16) can be substituted by infinity for all u, i. e., 



au 
PeO 

r c au 



» 1. 



As the u variable varies in the interval < u < 1, the slow cooling condition is 

r c a < p% < 1 . 

Then, with an exponentially small error, 



1 jl ( 1 

p res (u) - — {r c au) au r 

au \ au 



(r c au) a " T 1 + 



(5.17) 



(5.18) 



(5.19) 



This expression gives the probability that the spins in level u be in the up state at very low temperatures. In Fig. |^, 
the residual probability p ICS (u = 1) is plotted as a function of the cooling rate r c , for a = 1. The asymptotic result, 
given by Eq. (5. IE), is compared with the numerical integration of Eq. ( |5.16 ) withp e0 = 1/2, i. e., the system is taken 
initially at infinite temperature. The agreement is quite good up to r c ~ 0.1, which is not very small. 



P (u=l) 

* res^ 




FIG. 2. Residual probability for the slowest modes p re s(u = 1) as a funct ion o f the dimensionless cooling rate r c defined in 
the main text. The circles c orres pond to the numerical integration of Eq. (5.16), while the solid line is the prediction of the 
asymptotic calculation, Eq. (|5 . 1 9|) . A good agreement is observed up to r c ~ 0.1. 



The evolution of the average energy, for the cooling law (5.11), is given by 



e{t) - e e (T) = du [p(t; u) - p e 
Jo 



1 rPeo 
du / dp' e exp 



n 'au au 

"e ye 

r c au 



(5.20) 
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The residual energy can be easily computed by particularizing this expression for T — > 0, 



e rcs = / dup res (u) . 



(5.21) 



The integrand p ros (u), given by Eq. ( |5.1S| ), vanishes exponentially in the limit u — * 0. A standard Laplace analysis 
can be made, with the result 



Ercs ~ T 1 



The leading behavior is potential with r c , since 



1\ (r c a)~< 



In (r c a)~ 



(5.22) 



In e res ~ - In (r c a) 
a 



(5.23) 



which comes from the upper limit of integration, u = 1, corresponding to the largest relaxation time. The integral 
over the whole distribution function p les (u) gives a logarithmic correction |ln(r c et)|~ , which makes the residual energy 
smaller than the dominant term (rca) 1 / . This is due to the increasing behavior of p res (u) with u, p rcB (u) < p rcs (u = 1). 
In Fig. H the residual value of th e ene rgy i s plotted. The asymptotic expression, Eq. (5.22), is compared with the 
numerical results from Eqs. ( 5.21 ) and ( |5. 16 ). A good agreement is found up to r c = 0.01. It is worth noting that, for 
the cooling law considered, the logarithmic correction to the potential in r c behavior is not present for other simple 
models of structural glasses previously studied [^H-^5| . 




FIG. 3. Dimensionless residual energy e rcs as a function of the cooling rate r c . As in Fig. ^, the circles are from the 
numerical integration of Eq. (5.16), and the solid line is the prediction of the asymptotic analysis, given by Eq. (5.22). The 
agreement is good for r c < 0.01. 



A. Demarcation mode, Active temperature and glass transition 



The existence of non- vanishing residual properties is an indication of the departure of the system from equilibrium 
at low temperatures. Due to the divergence of the characteristic relaxation times r„, for low enough temperatures the 
system does not have enough time to relax towards equilibrium, and a kinetic phenomenon resembling the laboratory 
glass transition |22|-Q shows up. Next, we will try to understand the physical origin of this kinetic transition. 
Let us c onsi der again the time evolution of the probability distribution p(t; u) in a cooling process, as given by Eqs. 



( |5.5D and (|5.6| ). The integral in x(t,t';u), 



I{t,t';u) = / dt"a{T";u) (5.24) 
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is a measure of the average number of transitions occurring in level u in the time interval between t' and t. Conse- 
quently, a mathematical definition for the limit of slow cooling is that the condition 



l(t*,0;u) » 1 



(5.25) 



holds for all u, where t* is the ti me for which the temperature vanishes if extrapolated accordingly to the prescribed 
cooling law, i. e., T(t*) — 0. Eq. ( 5.25 ) guarantees that the system experiments a large number of transitions before 
getting eventua lly fr ozen, so that it has enough ti me to forget the details of the in itial condition. For the cooling law 
defined in Eq. ( 5.11 ) it is easily verified that Eq. ( [5.25 ) is equivalent to Eq. ( 5.17 ). 

If we are dealing with a slow cooling process, Eq. (5.25) implies that there is a time window over which 



l(t,0;u) » 1 , 



(5.26) 



This is the time regime we are int erest ed in. Let us analyze the behavior of x(t, t'; u) as a function of t' , to <t' < t, 
for a given time t such that Eq. (5.26) holds. The function I(t,t';u) changes from a very large value to zero when 
if goes from to t. Consequently, x\t,t';u) increases from practically zero to unity when t' moves in the above time 
interval. Let us define a time f/(t;u), prior to t, by 



I(t,tf; u) = 1 



(5.27) 



so that the average number of transitions taking place in level u in the time interval between t f and t equals unity. 
Then, x(t, tf]u) — e _1 , and the function x(t, t'; u) changes from zero to unity in a certain time interval around tf. In 
order to proceed, we will assume that this change takes place in the vicinity of tf very rapidly, as compared with the 
variation of the rest of the integrand of Eq. (|5.5|). Decomposing Eq. (5.5) in the form 



p(t;u) =p e {T) 



dt 



,d Pe (T)dT 



7 x(t,t';u) 



dT' dt 

the first integral is subdominant with respect to the second one and, moreover, 

.dp e (T')dV 



1 , ,d Pe {T')dT' . , . 

dt iJ^\ — >_ x (t,t';u) 

dT' dt' AV ' ' ' 



1 , ,dpJT) dT' . . 



dt' 



dT' dt 



dT' dt 

where Tf(t;u) is the temperature of the system at time tf(t;u), i. e., 

T f (t;u)=T[t f (t;u)}. 



7 =Pe{T) -p e [Tf(t;u)] , 



(5.28) 



(5.29) 



(5.30) 



Note that Tf(t;u) > T, since the time instant tf(t;u) < t. Substitution of Eq. (5.2£) into Eq. (5.28), and use of the 
above approximations yields 



p(t;u) =p e [T f (t;u)} 



The arguments leading from Eq. ( p.5[ ) to Eq. (5.31) are formally equivalent to assume that 

X(t,t';u) = e[t'-t f (t;u)], 



(5.31) 



(5.32) 



where Q(x) is Heaviside's step function. The result in Eq. (5.31) has a neat physical interpretation: the probability 
distribution of spins in level u at time t is given by th e equilibrium distribution corresponding to the temperature of the 
system at a prior time tf(t] u), defined in Eq. ( 5.27 ). The temperature Tf(t] u) is, therefore, the "Active" temperature 
of level u for the cooling program under consideration, in the sense used in the phcnomenological theories of glass- 
forming liquids 1 22 24|]. In the present context, each level has its own fictive temperature, so that there is not a 
unique equilibrium distribution describing the whole state of the system at a given temperature. As a consequence, 
the behavior of different macroscopic properties of the system may be quite different, depending on the modes which 
dominate for the evaluation of each property. 



The definition of the fictive temperature, Eq. (5.30), can be written as 

I{t*,t f \u) =I(t*,t;u)+I(t,t f ;u) = I(t*,t;u) + 1. 



(5.33) 



In the "high temperature" regime, where the number of tr ansitions between temperature T and T = in level u is 
very large, we can neglect the unity on the rhs of Eq. ( 5.33| ) and we get Tf(t; u) ~ T; the modes in level u remain in 
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equilibrium. On the other hand, in the "low temperature" limit, I(t*,t;u) becomes very small and Tf(t;u) tends to 
a limiting value, 



77 (tt) = hm T f (t; u) = T[t}(u)} , (5.34) 



where is the time for which 



I(t*,t*f,u)= dta(T;u) = l. (5.35) 



The low temperature region for level u corresponds to temperatures such that the average number of transitions 
remaining to spins in level u bef ore formally reaching T = is smaller than unity. Since the rate a(T; it) decreases 
as a function of it, see Eq. ( |3.14 ), Tj? (it) is an increasing function of u, i. e., the Active temperature T^(u) is larger 



the slower the modes. Then, the following picture of the kinetic glass-like transition appears. The relaxation modes 
associated to level u are in equilibrium for T > TJf(u), getting frozen in their equilibrium distribution corresponding 
to the Active temperature TJ(u) for T < Tj(u). A similar idea was first used in the context of the structural glass 
transition by Dyre p6| , who characterized the transition by means of the concept of demarcation mode. For a given 
value of the temperature T, the modes in levels with Tf(u) > T will be frozen, while those with TJ{u) < T will be still 
able to relax. The demarcation mode at temperature T, Ud(T), is defined by Tj(ud) — T. The modes with it > Ud(T) 
are frozen at this temperature, and modes with u < Ud{T) can still evolve. The glass transition can be understood 
as the movement of the demarcation mode from u = 1, corresponding to the largest relaxation time, to u = 0, which 
is the fastest mode. The temperature T g at which the glass transition begins for a given cooling program is that for 
which the slowest modes freeze, i. e., T g = T?(it = 1). 

The physical image developed in the above paragraph allows us to estimate the residual values in quite a simple 
way. For a given cooling law, first we calculate the limit value of the Active temperature T^{u). Then, the residual 
probability distribution is given by, 

Pres(lt) ~ p e [Tf(u)] , (5.36) 

and, consequently, the residual value of the energy is 

Ercs = / dup IeS (u) ~ / du p e [Tf (u)] . (5.37) 







As a tes t, we have considered the cooling process of the hierarchical model studied in this paper with the law given 



by Eq. ( 5.11 ). In the slow cooling limit, it is easy to show that 



T* f (u)^-^\ (5.38) 
1 \[n(r c a)\ 



i. e., the glass-like transition begins at 



2a 
In(r c a)| ' 



T 9 - mTTTTTF . ( 5 - 39 ) 



for which the slowest modes become frozen. In Fig. ^ we plot the evolution of the probability, as given by Eq. (5.13), 
of the slowest level u = 1, for a = 1 and r c = 10~ 4 . For these parameters, it is T g = Tj? (it = 1) ~ 0.217, which is seen 
to be a good value for the Active temperature. It is worth noting that the simple arguments used in this section lead 
to an estimate of the residual population 

. i 

p Tes (u) ~ p e [Tf(u)\ = (r c au) au , (5.40) 
that gives the same dependence of p Yca on r c as Eq. ( 5.19| ) , except for a factor of the order of unity. 
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FIG. 4. Probab ility for the slowest level, p(u = 1), as a function of the dimensionless temperature in a process wit h the 
cooling law of Eq. (5.11) and a cooling rate r c = 10 -4 (squares). The glass transition temperature T 3 , as given by Eq. (I5.39T ) 
is indicated. It gives a good estimate of the actual Active temperature for this level. The equilibrium probability p e , Eq. ( |2.5[ ), 
is also plotted (dotted line). 



In Fig. [5] the dimensionless average energy e, Eq. ( 5.20| ), is shown, for the same parameters as in Fig. [|. The 
predicted glass transition temperature T g ~ 0.217 is also indicated. It provides a good estimate for the beginning of 
the d epart ure fr om th e equilibrium curve. The residual value of the energy calculated from an asymptotic analysis of 
Eqs. $Mj and <EM) reads 



i _1 

, / \ (r c a) a 
dup ICS (u) ~ j 

In (r c a) a 



(5.41) 



Comparison of Eqs. ( 5.41 ) and Eq. ( [5.22] ) also indicates a good agreement, the difference being again a factor of the 
order of unity. We can introduce a global fictive temperature T j for the energy, as the temperature for which the 
equilibrium energy is the same as the energy of the system in the cooling process, when extrapolated to T — 0, i. e., 



.(57): 



(5.42) 



This global fictive temperature T ^ is related to the fictive temperatures of the levels Tj? (u) by 



e e (T f )= / dup e [T}{u)] 
Jo 



(5.43) 



The global fictive temperature is then a kind of average of the fictive temperatures Ti(u), but the probability 

distribution of the frozen state is not the equilibrium distribution at the temperature Tf. Each level freezes at 
its corresponding fictive temperature Tt(u), and the probability of the frozen state p re s(u) is approximately given by 

Eq. (5.40). For the values of the parameters considered in Fig. H, Eq . (5.41) leads to Ty ~ 0.175, which is also a good 
approximation for the value of the temperature at which Eq. ( 5.42 ) holds. Finally, let us stress that similar results 
for the fictive temperatures can be derived by means of the qualitative reasonings developed in Ref. p8|, and based 
on Hilbert's expansion. 



15 



10"' 




FIG. 5. Average dimensionless energy as a function of the dimensionless temperature T in a cooling process (squares). The 
cooling law and the cooling rate are the same as in Fig. ^. The dotted line is the equilibrium energy, Eq. (|2,6|). The predicted 
value for the glass transition temperature T g ~ 0.217 is indicated. Also, the global Active temperature Tf, calculated from Eq. 
( 5.41), is shown. It gives quite a good estimate of the temperature for which the equilibrium energy equals the residual energy. 



VI. HEATING PROCESSES 



This section will be devoted to the study of heating processes after the system has been previously cooled down 
to very low temperatures. Then, the system is not initially at equilibrium, but the initial probability distribution 
p(to;u) corresponds to the final state reached in the cooling process. It is convenient t o in troduce a time th < to such 
that T(th) = if the heating law is extrapolated to times shorter than to- From Eq. (4.2) it follows that 



p(t;u) 



\p(t ;u) -p e (T )]x(t,t ;u) 

, .dp e (T') dT' , , 
dt , -^ r ^— X (t,t';i 



Pe(T) 



dT' dt' 



(6.1) 



where x(ti,t2]u) was define d in Eq. (5.6), and T = T(to) is the initial temperature of the heating process. The 
first term on the rhs of Eq. ( |3.l| ) represents the decay of the initial nonequilibrium condition. For T(t) in the high 
temperature region it is x(t,to;u) <C 1, and this contribution can be neglected. This implies that p{t;u) reaches a 
behavior which is independent of the initial condition or, equivalently, independent of the previous cooling program. 
Moreover, the lower limit to in the integral can be replaced by th, with a relative error that decreases as t increases. 
Therefore, p{t; u) tends to the "normal" solution 



p N {T;u) =p e (T) 



dt , , J —x(t,t';u). 
dT' dt' AV ' ' ' 



(6.2) 



Any arbitrary pit; u) approaches pn(T; u) for long enough times, corresponding to high temperatures. The existence 
of the normal solution for heating processes is not restricted to the hierarchical models considered in this paper, 
but it has been established for a quite general class of systems, whose dynamics is described in terms of a master 
equation |27j]. In this context, the existence of the normal solution in our simplified, mean-field description of the 
hierarchical models, provides a consistency test of the approximations we have introduced starting from the master 
equation formulation. However, it cannot be assured that the special choice of the initial conditions leading to the 
normal solution in our mean field type approximation remains the same for a more exact analysis. 

In th e limi t of high temperatures, pn(T; u) must be closely related to Hilbert's expansion solution pu {T; u), analyzed 
in Sec. IV A. In fact, since pn{T;u) does not refer to any particular initial condition, it is to be expected that the 



normal solution coincides with Hilbert's expansion in their common range of validity. However, the normal solution is 
not restricted to near equilibrium situations, and the system can approach the normal curve in a temperature range 
for which Hilbert's solution is not accurate. This point will be clearly illustrated below. 
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For the sake of simplicity, we will consider that the system is submitted to a cooling process given by Eq. (5.11) 
followed by a heating process of the form 



dp e 
dt 



where > is the heating rate. Using this heating law, we can express Eq. (6.1) in the simpler way 



p au _ p au \ CP" / p au _ p'au 

p(t; u) = \po(u) - p e0 ] exp ( e e0 ) +p e - \ dp' e exp 



r h au 



Here we employ the same notation as in Eq. (5.12). Similarly. 



Pn{T; u)=p e - j dp' e exp - 

Jo V r hau 



au _ lau 
-Pe re 



The high temperature regime for level u is given by 



V" 



> 1 



Thau 



and, in this limit, Eq. (|6.5[) reduces to 



p N (T;u) ~p e - r h p i 



1 — au 



(6.3) 



(6.4) 



(6.5) 



(6.6) 



(6.7) 



that agrees with Hubert's solution. Eq. ( 4.11 ), particularized f or th e heating law ( |6.3| ). 

Figure [3] shows the evolution of p(t; u = 1), as given by Eq. (6.4), in a heating process with a 



1 and rh = 10 



Two different initial conditions have been considered, corresponding to slow previous coolings of the system with 
r c = 10~ 4 and r c = 10~ 6 , respectively. In the figure, it is seen that both heating curves approach the normal 
solution, reaching it in a region where the normal curve represents a clear nonequilibrium state. Hilbert's expansion 
approximation, Eq. (6/7), is also plotted. It provides a good description for the linear correction around equilibrium 
of th e normal curve, but fails for low enough temperatures. In fact, for the values of the parameters considered, Eq. 
(3.7) becomes negative for T < 0.217. Therefore, in general, Hilbert's solution cannot be used to estimate the normal 
curve over the whole range in which the system is well described by the normal solution, which depends on the details 
of the previous cooling process. In the figure, for the smallest cooling rate r c = 10 -6 , the system reaches the normal 
curve for a temperature at which Hilbert's expansion is not valid. This illustrates how the normal solution is relevant 
for far from equilibrium states, in which a linear theory in the deviations from equilibrium is not accurate. 
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FIG. 6. Evolution of p(u = 1) in a heating process. The dotted line is the equilibrium curve. The heating law is Eq. 
(|3.3|), with a heating rate th = 10~ 4 . The normal curve has been evaluated both by using Eq. (6.5) (dashed line) and 
Hilbert's expansion, Eq. (6.7) (solid line). Two different initial conditions for the heating process have been considered, both 
corresponding to previous coolings down to very low t empe ratures but with different rates, namely r c — 10 -4 (squares) and 



10 (circles). The cooling law was given by Eq. (5.11) in both cases 
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In the experiments, either with supercooled liquids or with granular materials, the distribution function cannot be 
directly measured, but instead the average values of the relevant physical properties. In our model, we can consider 
the dimensionless energy e(t) |50|, whic h is calculated by integrating p(t;u) over u. For the specific heating law we 
are considering, it follows from Eq. (3.4) that 



e(t) 



du [ Po (u) - Pe0 ] exp ( -£s E-o ) + £e{T) 



du / dp' e exp 

n Jp e0 



r h au 

n 'au 



r h au 



(6.8) 



There is also a normal curve for the energy e(T), that can be obtained by integration of the normal probability 
distribution, Eq. (6.5), i. e., 



du / dp' e exp 



Pe 



r h au 



(6.9) 



For long enough times, independently of the initial probability distribution po(u), e(t) will approach ejv(T), as a 
consequence of the tendency of p(t;u) towards Pjy(T ; u). The time regime in which e(t) practically agrees with 
£jv(T) corresponds to the condition given by Eq. (6J3) being verified for all u, i. e., when it is fulfilled by the slowest 
modes, u = 1. Therefore, the appro ach to the normal curve is controlled by the slowest l eve l u — 1. For very high 
temperatures, Hilbert's result (4.13) holds and, particularizing for the heating law of Eq. (6.3), we get 



e N (T) ~ e H (T) = p e - r h p ( 



Pe ~ 
a\ lnp e 



(6.10) 




T 

FIG . 7. Hysteresis of the dimensionless energy e in a thermal cycle. As in the previous figures, the cooling law is Eq. 
(5.11) and the heating law is Eq. (|6.S|), The cooling and heating rates are the same, r c — rn = 1CP 4 . The dotted line is the 
equilibrium energy, and the dashed line is the normal curve for the heating process. The arrows over the solid lines indicate 
the variation of the temperature in each process. In the heating process, the system approaches the normal solution, crossing 
the equilibrium curve. 



An important feature of the normal curve is that 

e N (T) < e e {T) , 



(6.11) 



for all temperatures. This property, together with the inequality Eq. (5.4) for cooling processes, explains the hysteresis 
effects exhibited by the model when submitted to a thermal cycle. These effects are similar to the experimentally 
observed behavior in glass-forming liquids p2[ and other complex systems, as granular materials pq |. In Fig. |?] a 
hysteresis cycle of the energy is shown for r c = = 1CP 4 . In the cooling process the energy is greater than the 
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corresponding equilibriu m valu e, while in the heating process the tendency of the system to approach the normal 
curve, which verifies Eq. ( |6.11 ), makes the system overtake the equilibrium curve. Only for very high temperatures, 
where Hilbert's expansion is accurate, the normal solution tends to the equilibrium one from below. Comparison of 
Figs. | and | shows that the separation of the normal solution for the energy from the equilibrium curve is smaller 
than that of the slowest level probability p(u = 1). This is due to the contribution of the other modes, which reach 
their equilibrium values for lower temperatures. 

In our opinion, the explanation of the hysteretic behavior observed in real systems must be similar to the above 
discussion. This idea is supported by the fact that in other models based on a master equation formalism, analogous 
results have been found |2^,|3l],|35| . Also, Hilbert's expansion, although of limited validity, leads to the conclusion that 
cooling and heating curves are at opposite sides of the equilibrium curve. 



VII. DISCUSSION 



In this paper, we have analyzed a simple spin model with hierarchically constrained dynamics. The spins are 
classified into levels, and a spin in level n + 1 is able to relax only if a certain number of spins in level n are in the 
up (excited) state. Starting from a master equation formulation of the dynamics, a mean-field approximation was 
introduced. The result is a generalization of Palmer's et al. model Q, with a closed equation for the evolution of the 
probability of finding a given spin in the up state. Each level has a characteristic relaxation rate, which is a function 
of the temperature of the system. Two sets of parameters characterize the model: the number of spins in level n, N n , 
and the number fi n of spins in level n involved in the facilitation of the relaxation of a spin in level n + 1. We have 
chosen the simplest possibility, namely that \i n and N n are proportional to each other, which leads to a relaxation 
behavior that is independent of the level populations N n . In relaxation processes at constant temperature, the system 
displays linear logarithmic decay jllj . This is a characteristic feature of the behavior of a wide variety of complex 
systems, including structural glasses, spin-glasses, protein models and powders. 

A key point in our approach to the problem of hierarchically constrained dynamics is the introduction of the 
"mean-field" approximation, which allows us to reduce the initial, rather involved, problem to a solvable one. This 
approximation is based on the physical idea that hierarchical constraints render the characteristic time scales of the 
different levels clearly separated. This is indeed the case in the low temperature region, since the ratio t„/t„_i of the 
characteristic times of levels n and n — 1, as given by Eq. (3.5), diverges for T — > 0. 

The model has also been used to study processes in which the temperature changes in time in an arbitrary way. 
The general solution for the time dependent distribution function can be explicitly written, and the average energy 
has a form resembling that of Narayanaswami's phenomenological theory of glass- forming liquids |25|. By means of 
Hilbert's expansion we have constructed an approximate solution for the probability distribution. This expression is 
valid in the limit of very high temperatures, and for situations where the system is in the linear around equilibrium 
region. The behavior of Hilbert's solution is formally identical to the one previously found for a very general class 
of systems described by master equations p8fl , despite the mean-field character of the simplified model considered 
here. Hilbert's expansion predicts the existence of hysteresis effects when the system is first cooled down to low 
temperatures and, afterwards, reheated to high temperatures. This is because the average energy is at opposite sides 
of the equilibrium curve for cooling and heating processes. 

Another point we have addressed is the evolution of the system in continuous cooling processes. A phenomenon 
similar to the laboratory glass transition shows up, and the system departs from equilibrium at low temperatures. The 
magnitude of this separation can be measured by the values of the residual properties, which have been analytically 
computed. A simple but physically appealing argument is presented, in order to understand the origin of this glassy 
behavior. Each level in the system becomes frozen at the equilibrium value corresponding to a temperature, called the 
fictive temperature of the level, such that the average number of transitions per spin in that level until reaching T = 0, 
following the prescribed cooling program, equals unity. A similar argument has been previously used in other models 
p8| , pl|j35| . Here, an analytical derivation is presented. This provides a theoretical basis for the concept of fictive 
temperature, and clarifies the accuracy of the results following from the qualitative argument. We have compared 
them with the values of the residual properties obtained numerically as well as with those following from asymptotic 
analysis calculations. 

Finally, heating processes have also been analyzed. More specifically, we have considered heating processes following 
a continuous cooling of the system down to very low temperatures. In the description of heating processes, the so- 
called normal curve plays a fundamental role. The hysteresis effects observed when the system is first cooled and 
afterwards reheated appear because of the trend of the system to approach the normal solution, along the heating 
evolution. This is analogous to previous results found in models for structural glasses |3l],|35| and for granular systems 
p8[. Thus, we think it would be worth investigating whether a similar curve does exist for real complex systems 
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showing this kind of hysteretic behavior. 

Although a particular simple model with hierarchically constrained dynamics has been considered in this paper, 
most of the physical ideas developed seem very general. The validity of Hilbert's expansion for high temperatures, the 
natural appearance of the concept of Active temperature, leading a physically appealing description of the laboratory 
glass transition, the existence of the normal solution and its fundamental role in explaining the hysteresis effects, are 
results that are not restricted to the present model. In fact, similar results appear in a quite general class of systems 
described by master equations Of course, the details of the dynamical behavior depend on the specific 

model we are dealing with. On the other hand, the general picture developed here, where the heating and cooling 
experiments in glasses appear as purely kinetic and relatively simple to understand, might not be valid for real and 
much more complex systems. 
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